% 增长区域的代码
function[c,f,s]=pdefun(x,t,u,ux)

c = [1;1];

% 固定 Exm1
% ρ(t)=1, d1=0.35，d2=0.25, theta=2.4, b=0.21, c=0.4, m=0.5, r=1.7, f=0.3
% y(1) = 2.4 * u(1) - 0.21 *(u(1))^2 - 0.4 * u(1) * u(2) / (0.5+u(1));
% 
% y(2) = 1.7 * u(2) - 0.3 * (u(2))^2 / (0.5 + u(1));
% F=[y(1);y(2)];
% 
% s=F.*[1;1];
% 
% f = [(0.35*ux(1));(0.25*ux(2))];


% 固定 Exm2
%ρ(t)=1, d1=0.009，d2=0.2, theta=2.4, b=0.21, c=0.4, m=0.5, r=1.7, f=0.3
% y(1) = 2.4 * u(1) - 0.21 *(u(1))^2 - 0.4 * u(1) * u(2) / (0.5+u(1));
% 
% y(2) = 1.7 * u(2) - 0.3 * (u(2))^2 / (0.5 + u(1));
% 
% F=[y(1);y(2)];
% 
% s=F.*[1;1];
% 
% f = [(0.009*ux(1));(0.2*ux(2))];


% 增长 Exm1
% ρ(t)=(exp(5*t))/(1+0.2*(exp(5*t)-1)), d1=0.35，d2=0.25, theta=2.4, b=0.21, c=0.4, m=0.5, r=1.7, f=0.3

% p = (exp(5*t))/(1+0.2*(exp(5*t)-1))
% 
% y(1) = 2.4 * u(1) - 0.21 *(u(1))^2 - 0.4 * u(1) * u(2) / (0.5+u(1));
% 
% y(2) = 1.7 * u(2) - 0.3 * (u(2))^2 / (0.5 + u(1));
% 
% 
% F=[y(1);y(2)];
% 
% s=F.*[1;1];
% 
% f = [(0.35*ux(1))/(p^2);(0.25*ux(2))/(p^2)];




% 增长 Exm2
% ρ(t)=(exp(5*t))/(1+0.2*(exp(5*t)-1)), d1=0.009，d2=0.2, theta=2.4, b=0.21, c=0.4, m=0.5, r=1.7, f=0.3

p = (exp(5*t))/(1+0.2*(exp(5*t)-1))

y(1) = 2.4 * u(1) - 0.21 *(u(1))^2 - 0.4 * u(1) * u(2) / (0.5+u(1));

y(2) = 1.7 * u(2) - 0.3 * (u(2))^2 / (0.5 + u(1));


F=[y(1);y(2)];

s=F.*[1;1];

f = [(0.009*ux(1))/(p^2);(0.2*ux(2))/(p^2)];


end






% y(1)=u(1)*(1.8+0.3*abs(sin(pi*t))-(0.3*abs(sin(pi*t)))*u(2)*inv(1+0.4*abs(sin(pi*t))*u(1))-u(1));
% f = [0.14*inv(p^2)*ux(1);0.10*inv(p^2)*ux(2)];
% y(2)=u(2)*(1+0.5*abs(sin(pi*t)))+(0.2*abs(sin(pi*t))*u(1)*inv(1+0.4*abs(sin(pi*t))*u(1))-u(2));